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Abstract 

We introduce a method to obtain self-consistent, axially symmetric, 
thin disklike stellar models in the first post-Newtonian (1PN) approx- 
imation. The models obtained are fully analytical and corresponds to 
the post-Newtonian generalizations of classical ones. By introducing 
in the field equations provided by the 1PN approximation a known 
distribution function (DF) corresponding to a Newtonian model, two 
fundamental equations determining the 1PN corrections are obtained, 
which are solved using the Hunter method. The rotation curves of 
the lPN-corrected models differs from the classical ones and, for the 
generalized Kalnajs discs, the 1PN corrections are clearly appreciable 
with values of the mass and radius of a typical galaxy. On the other 
hand, the relativistic mass correction can be ignored for all models. 

PACS numbers: 04.25.Nx, 98.10,+z 

1 Introduction 

The stars observed in the universe tend to cluster in huge self-gravitating 
systems, being the galaxies among the most noticeable and studied of them. 
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Galaxies can be described by models whose mass is distributed symmetri- 
cally about an axis to a finite distance of it i.e., the galaxy radius. Par- 
ticularly useful models are those restricted to a plane or thin disk, since 
for many galaxies like the Milky Way its height is small compared to its 
radius. For galaxies in general, the average time between collisions or in- 
dividual meetings (mean collision time) is greater than the system's life 
time. Many of the current models take Newton's law as their field equations 

mmOHISElElElUnnilll]. However ' many models have been developed 
recently under general relativity p21 \T3\ EH E3 EH El [THJ [H] , being one 

of the main motivations for including corrections made by general relativity 

the actual incompatibility between the rotation curves of theoretical models 

and the ones observed. 

The post-Newtonian approximation introduces general relativity through 

a series over the speed (v/c). For Minkowski's background metric, the first 

post-Newtonian corrections (1PN) are included taking the terms 

2 4 

g 00 ~ -1+ 9oo + 9oo, (1) 
9i0 ~ 9i0, (2) 

2 

9ij ~ 9ij, (3) 

where the upper index denotes the order of the power of (v/c) of the term. 
Using harmonic coordinates, the 1PN order potentials are defined as 

goo = -2<P/c 2 , (4) 

9i0 = d/c 3 , (5) 

500 = -2(0 2 + VO/c 4 . (6) 

Henceforth we consider only the stationary case, so the explicit time-dependent 
terms along with the potential vector £ disappear in all equations. 
The stationary field equations in the 1PN approximation are 



4ttG 



o 



VV = —rf- T 00 , (7) 

2 2 

V 2 V> = 4^G(T 00 + T u ), (8) 

where, in the classical limit, potential ip vanishes and eft tends to the New- 
tonian potential cpN- Whereas that 

^ = -V0 - ^ [V (2(f) 2 + ij>) + 4v(v • V0) - v 2 X7<p] (9) 
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is the stationary equation of motion. 

A complete statistical description is achieved by knowing the distribution 
function (DF) of the system. The DF satisfies a continuity equation in the 
phase space called the Boltzmann equation. In the 1PN approximation, 
the collisionless Boltzmann equation (CBE) for a many identical particles 
system is given by [20] 



dF d<f>dF 1 (dcj), Al 2 , d(p , , dtp \ dF 



v % — r — r — (4</> + V s ) — uV H -7 r = 0. (10) 

dx l dx 1 dv l c 2 \ dx 1 dxi dx % J dv l 

According to Jeans theorem [21J, the solution of colissionless Boltzmann 
equation is any function of the integrals of motion. Now, it is easy to verify 
that two isolated integrals a 1PN system with axial symmetry are given by 

E=~v 2 + <!>, (11) 

where 

26 2 + ib 

<S> = 4> + v 2 - , 12 

and 

L z = Rv^e-^ c2 (13) 

which can be interpreted as the 1PN generalizations of energy and the z 
component of angular momentum, respectively. 

In addition, the DF must satisfy the condition of self-consistently gener- 
ating the macroscopic mean values. In the post-Newtonian approximation, 
the following components of the energy-momentum tensor are need, 

o 



x,t) = c 2 J F(x,v,t)d 3 v, (14) 
T 00 + T li = 2 j (t> 2 (x, t) + <Xx, t))F(x, v, t)d 3 v. (15) 



Therefore, in the post-Newtonian approximation, self-consistent equilibrium 
models are defined by two scalar potentials, eft and ip, together with a DF that 
satisfies 1PN CBE and relations (|14p - (I15p generating self-consistenly the 
1PN components of the energy-momentum tensor. In this paper we present 
a method to obtain post-Newtonian axially symmetric equilibrium models. 
The method uses thin disk models, allowing to solve the two differential 
equations with the Hunter's method [22], as shown in the next section. The 
solution of the field equations is obtained considering equations in vacuum, 
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in which case, the energy-momentum tensor vanishes and the content of 
matter is expressed as a boundary condition on the fields in the disk, as 
shown below. Finally, in section [3] the first axially symmetric models in the 
1PN approximation are presented. 

2 Formulation of the method 

For thin disks of finite radius the components of the energy-momentum 
tensor can be written as 



T 00 = c 2 Z(R)5(z), (16) 

2 2 

r oo + T u = a ( R )S(z), (17) 



for < R < a, where 5 is the Dirac delta function, and being zero for R > a. 
Therefore, the field equations reduce to two Laplace equations for the fields 
(j) and ip. It is demanded that the fields are even functions of z 

(t>(R,z) = <j>(R,-z), ip(R,z) = il>(R,-z), (18) 

and therefore, that the first derivatives with respect to z are odd functions 
of z. 

Using Gauss' theorem with (|16p and (117j) gives, 



=0+ 



2vrG V dz 



The problem is defined with the following boundary conditions: fields vanish 
at infinity, and at z = its derivatives depend on the energy-momentum 
tensor components according to (|19|) and (|20p for < R < a and vanishing 
for R > a. Applying Hunter's method to each of the 1PN field equations, 
one can obtain exact analytical expressions for the potentials. The Hunter's 
method consists on obtaining solutions of the Laplace equation in oblate 
spheroidal coordinates. The oblate coordinates are related to the cylindrical 
by 



r = ava+aa-r? 2 ), (21) 

z = a^r], (22) 
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where < £ < oo and —1 < rj < 1. The disk is placed at £ = 0, where, 
rf = \- R 2 /a 2 . 

Following Hunter [22], the general solution of each Laplace equation 
satisfying previous boundary conditions can be written as 

oo 

ri) = ~J2 A 2nq2n(0P2n(v), (23) 
n=0 

for the (j) potential, and 

oo 

W,V) = ~J2 B ^<l2n(OP2n(v), (24) 

n=0 

for the i/j potential, where Am and B 2 n are the constants required for each 
model, P2n{ji) are the Legendre polynomials and q2n(C) = i 2n+1 Q2n{iQ are 
Legendre functions of second kind. Note that when using classical models, 
the expression for (j) can be written taking constants of the form 

A 2n = C 2n + D 2n /c 2 ) (25) 

where the C2 n constants define the Newtonian potential <Pn and the con- 
stants Z?2n define the correction <Ppn- So that taking the limit c — » oo, 
<P = 4>N + <t>PN is reduced to the Newtonian part only. The corresponding 
expressions for E and a in oblate coordinates are 



£ = j^ Mn(2n + ^n+lW^nM, (26) 

1 x 

a = r— ^- Y\ B 2n(2n + l)g 2 n+l(0)P 2 n(r?). (27) 
2-iraGri ^-^ 

n=0 

Accordingly, E can also be writen as the sum of a Newtonian part and a 
post-Newtonian correction: E = Ejv + Epjv- 

In order to obtain self-consistent models of axially symmetric thin disks 
in equilibrium, it is used a DF of the form F = f(R, v R , v v )S(z)5(v z ), which 
is zero for R > a. The DF reproduces E and a through the equations 



E(i?) = J J f{R,v R , Vip )dv R dv v , (28) 
a{R) = 4j J Ef(R,v R , Vcp )dv R dv v -2ct> N Z N , (29) 
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where E is the integral of motion defined by (jlip . As in the Newtonian case, 
one could define a relative energy and a relative potential as 



e = -E + $o, 
* = -$ + $ , 



(30) 
(31) 



where <l?o is a constant that is chosen so that e and are always positive, 
i.e., such that / > for < e < ^. 

The method to obtain self-consistent post-Newtonian models is to take 
(|28p with (|26p as the first fundamental equation and (|29p with ()27p as the 
second fundamental equation. In accordance with the order of magnitude of 
the approach, the second fundamental equation is solved using the Newto- 
nian terms and then the first fundamental equation is solved using the terms 
up to 1PN order. Finally, note that for the second fundamental equation, 
the constants can be obtained directly by 



where the integral depends on the particular model. On the other hand, in 
the first fundamental equation, the constants also appear on the right side 
of the equation (hence, it is necessary to obtain the constants differently for 
each model) and, in general, we can not use an explicit expression as done 
above. However, it is found that for all models treated it is finally possible 
to obtain expressions analog to (|32p for the first fundamental equation, see 
(|57j) - (JSnj) and (JHD - ff5]) . The DF in the 1PN approximation presents the 
same functional dependence on the integrals of motion that in the Newto- 
nian models, that is, one uses the same DF but now with the 1PN energy 
and angular momentum. Of course, this is for the sake of correspondence 
with the Newtonian limit. Note that the correspondence principle must be 
satisfied by both the DF and its integrals. 

It is possible to obtain alternative expressions for S and a considering in 
particular the case where the DF depends on a linear combination of energy 
and angular momentum called Jacobi's integral, J = £IL Z — E. Jacobi's 
integral is interpreted classically as the energy measured from a frame of 
reference rotating with angular speed [21]. In terms of the relative energy 
([30]) . we can write Jacobi's integral as 



where ^ e (0) is the 1PN relative-effective potential evaluated in r] = 0, de- 
fined by 




(32) 



J = e + fiL z -¥ e (0), 



(33) 




(34) 
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The Jacobi's integral takes values between zero and J m ax, with 



1 



Jmax = V - ^ 2 a 2 ri 2 11 2 
Now, given that 2nd J = dvRdv^, the relation (f28l) can be written as 

pJmax 

Z = 2tt f(J)dJ. 
Jo 

Then, by using the expression 
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v^fdvRdv^ = 2-k(v v ) 



Jmax 



f(J)dJ, 



(35) 



(36) 



(37) 



the second fundamental equation can be written as 



r Jmax 

a = 2(2<J» - V 2 a 2 - <p N + 2afi\A - »7 2 (*V>) S JV - 8tt / Jf(J)dJ, (38) 



with S^v given by (|26|) with Am = C2n- 

Finally, the circular speed necessary for the rotation curve can be ob- 
tained considering ([9]) for circular orbits. In that case, the circular speed is 
equal to and perpendicular to the gradient of the field 4>. Therefore, the 
1PN equation of motion ([9]) reduces to 



R 



c 2 3R 



z=0 



d_ 

dR 



+ 



2<p 2 + V 



(39) 



,2=0 



Then, in accordance with the 1PN order of approximation, the expression 
for the circular speed is 



IR 



dR 



1 + 



4(p R 



OR 



+ 



Rdip 
l?d~R 



(40) 



2 = 



Note that in the limit c — ► oo, the previous expression reduces to that of 
Newtonian theory v v = \J Rd<f>N / dR. 



3 Application to the Generalized Kalnajs Disks 

We will now apply the previous formalism to the family of Generalized 
Kalnajs Disks (GKD) introduced by Gonzalez and Reina in [9J. This family 
is characterized by mass densities of the form 

( m ) _ 3M 2m-l (AA \ 

Ljv '2^ ' (41) 
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where the index m of the model is any positive integer. The potential cj)ffi 
is given by (|23|) by taking the Newtonian limit in ([25p and using 

c [m) = MG^/ 2 (4n + l)(2m + l)! 

2n a2 2m + 1 (2n + l)(m-n)!r(m + n + |)g 2 n+i(0) 

for n < m, and = for n > m. 

3.1 The 1PN model for the m = 1 GKD 

The first disk of the family, the disk with m = 1, is the well known Kalnajs 
disk [6], with the mass distribution 

^ N ~ 2^ (43) 

Then, as was shown by Kalnajs [23], the DF depends on the Jacobi's integral 
as 

/(J) = w[ 2( ^-^ 2)j r 1/2 ( 44 ) 

where VL^ = 3irGM/4a 3 . Now, for the Kalnajs disk it holds that (v^) = UR, 
so that the disk spins like a rigid body. 

Inserting expressions (f4*3j) and (|4"4"|) into the second fundamental equation 
(I38p . and integrating the DF, one easily obtains the following system of 
equations 

B -B 2 + B 4 = 0, (45) 

27irC 2 M 2 

3B 2 -10B A = QMGati 2 — 5 , (46) 

4a z 

Bi = ^G 2 M 2 _ 3QGMn 2 a (47) 
140or 

Likewise, for the first fundamental equation we have 
^ 9ttG 2 M 2 Q 2 

D »- D2 + Dl " to>m-a>y <48) 

15 

jD 2 + #D 4 = B 2 -—B 4 , (49) 

o 

X D 4 = |S 4 - - GMtfa, (50) 



S 



where 

24a 3 (0 2 _ O 2 ) - 9itMG , , 

7 = ' (51) 

135ttMG/8 - 80a 3 (Og - ^ 2 ) 

* " MG ' (52j 

70[128a 3 (fi 2 - fi 2 ) - 27ttMG] 

X ~ SMttMG ' (53) 

Therefore, we have a system of linear equations with an upper triangular 
matrix for the constants of each potential, tp and 4>pn- 
Solving for the constants explicitly, we obtain 

B = 2aGM(TT-l)n 2 -^^I^ll, (54) 

B 2 = faGM(7vr - 10)O 2 - 3G2M ^g l7r ~ 2) , (55) 

d _ WGM ( GMir 5a 2 n 2 \ (rc\ 

^ 4 - ~l5a~ \~8^ S - J' 







"8a 2 (a 3 ft 2 -3GA/7r) 448a 2 (128a 3 H 2 -357GA/tt) 

. 3G 2 7r(l+807r)Af 2 _ 9(28G 3 A/ 3 7r 3 -95G 3 Af 3 tt 2 ) 
320a 2 28a 2 (8a 3 f2 2 -21GAf7r) ' 



(57) 



D _ _ 3G 2 7r(77r-13)M 2 _ 9(28G 3 Af 3 7r 3 -95G 3 M 3 7r 2 ) 

2 28^ 28a 2 (8a 3 n 2 -21GAf7r) ' ^ ' 

n _ 27G 2 M 2 7r(l84C 2 a 3 +39GM7r) (c .q\ 
4 ~~ 140a 2 (128a 3 Q 2 -357GA/7r) ' { ' 

which defines the potential ip, and the correction 0pjv, respectively. Hence 
we can get all the parameters of interest, in particular, we observed that the 
mass correction, Spjy, is negligible compared with the classical mass S^r. 
In contrast, the 1PN rotation curve obtained with (|40p . is visibly separated 
from the Newtonian one from a certain radius, the difference being maximum 
at the edge of the disk (see Fig. [1]), where the difference between both reaches 
10.3% approximately. Those have been obtained with the typical values of 
a galaxy such as the Milky Way. 

3.2 The 1PN model for the m = 2 GKD 

The procedure for the second disk is pretty much the same that for the first 
one. Again, the DF has a simple form when written as a function of Jacobi's 
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Figure 1: First we plot the Newtonian and 1PN mass density for the m = 
1 GKD. The two curves seem to be superposed meaning that the mass 
correction is negligible. Then we plot the Newtonian (dashed line) and 1PN 
(full line) rotation curves for the same disk. The 1PN corrections are clearly 
notorious reaching a maximum at the border of the disc. Parameter values 
are: M = 4 x 10 40 %, a = 3 x 10 20 m, J] = 2x lO~ 13 Hz. 



integral, 



2M ( 10a 3 \ 1/4 



for a mean circular speed (v v ) = Ml = yJlbnGM /32a 3 R. This DF could be 
easily integrated to self-consistently obtain the mass density of the model, 
which also can be obtained from (1411) with m = 2, 



(2) 5M / i? 2 \ 3/2 5M7? 3 

sV= wr7j =i^p- (61) 

the associated gravitational potential is given by the Newtonian limit 4>n of 
(1231) with (@2D and m = 2. 

Again, we replace the Newtonian terms into the second fundamental 
equation (f38|) to obtain: 

4 4 

^ 5 2n (2n + l)g 2n+1 (0)P 2n (r ? ) = C 2i?? 2t , (62) 

n=0 n=0 
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where Co = C2 = and 

C 4 = (25^ 2 10!G 2 M 2 )/(160a 2 ), (63) 



C 6 = -(75ir^W\G 2 M z )/{160a z ) (64) 

C 8 = (75^ 2 10!G 2 M 2 )/(2240a 2 ). (65) 

Multiplying (|62|) with a Legendre polynomial, integrating with respect to r], 
and using the orthogonality properties of the Legendre polynomials we get 

B _y ^C 2t 2- 2 ^(4n + l)r(2n + l) 

^ ^<Z2n+i(0)(2n + l)r(i-n + l)r(i + n + 3/2)' 1 ' 

We also could have used the expression ([32]) , 

After rearranging the terms, the first fundamental equation can be writ- 
ten as 

4 

Y J { D ^2nP2n(v)+q2n(0)P2n(0)}-B 2n q 2n (0)P 2n ( V )-C 2nV 2n } = 0, (67) 
n=0 

where, 

??2n = vr(2j + l)/(32a 2 (7 2n+1 (0)-g 2n (0)), (68) 

C = (675^ 2 G 2 M 2 )/(4096a 2 ) +^(0,0), (69) 

C 2 = -(1575vr 2 G 2 M 2 )/(4096a 2 ), (70) 

C 4 = -(1125vr 2 G 2 M 2 )/(2048a 2 ), (71) 

C 6 = -(2025vr 2 G 2 M 2 )/(4096a 2 ), (72) 

C s = (2025vr 2 G 2 M 2 )/(8192a 2 ). (73) 

Integrating and using orthogonality relations, the constants D 2n are given 
by 

n = A v ^C 2 ,2- 2i - 1 (4n + l)r(2j + l) B 2n q 2n (0) 
2n ^2n(2n + l)T(i -n + l)T(i + n + 3/2) + tf 2n ' 1 J 
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Figure 2: We plot the Newtonian and 1PN mass density, as well as the 
Newtonian (dashed line) and 1PN (full line) rotation curves for the m = 2 
GKD, with parameter values: M = 4 x 10 40 %, a = 3 x 10 20 m, 0, = 2 x 
lQ~ l ^Hz. As we can see, the two graphics of mass density appear to be 
superposed, while that the 1PN rotation curve behaves completely different 
to the Newtonian one. 



for n > 0, and 



Dn 



#0 + tt/2 



fTiTlWTsTol -2^^^(0)^2,(0) 



4 

§~r(i + i)r(i + 3/2) 



(75) 



which defines the correction 4>pn to the Newtonian potential. The mass 
densities and rotation curves are plotted at Figure [51 note that the rota- 
tion curve is not only cuantitative, but cualitatively different, presenting its 
maximum value more closely to the center of the disk and going upwards at 
the border. This is due to the fact that the 1PN curve includes more terms 
than the Newtonian one. 
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